#ifndef HUNGRAY_H
#define HUNGRAY_H

#include<iostream>
#include<algorithm>
#include<fstream>
#include<stdio.h>

#define SIZE 10
using namespace std;

class HungrayAlgorithmBase
{
public:
    HungrayAlgorithmBase();
    virtual ~HungrayAlgorithmBase();
    virtual void runMunkres() = 0;
};

class HungrayOrigin:public HungrayAlgorithmBase
{
public:
    HungrayOrigin()
    {}
    ~HungrayOrigin(){}
    void runMunkres()
    {
        bool done = false;
        while(!done)
        {
            showCostMatrix();
            showMasjkMatrix();
            switch(step)
            {
            case 1:
                step_one(step);
                break;
            case 2:
                step_two(step);
                break;
            case 3:
                step_three(step);
                break;
            case 4:
                step_four(step);
                break;
            case 5:
                step_five(step);
                break;
            case 6:
                step_six(step);
                break;
            case 7:
                step_seven(step);
                done = true;
                break;

            }

        }
    }
public:
    int C[SIZE][SIZE];// cost matrix [row] [column]
    int RowCover[SIZE];
    int ColCover[SIZE];
    int M[SIZE][SIZE];//  mask matrix
    int path[SIZE][SIZE];
private:
    /* variables */
    int step;
    int nrow, ncol;// row and column of the cost matrix
    int path_row_0;
    int path_col_0;
    int path_count;
    int MaxValue;

    void showCostMatrix();
    void showMasjkMatrix();
    void step_one(int step)
    {
        int min_in_row;
        for (int r =0;r < nrow; r++)
        {
            min_in_row = C[r][0];
            for (int c = 0; c < ncol; c++)
            {
                if(C[r][c] < min_in_row)
                {
                    min_in_row = C[r][c];
                }
            }
            for (int c = 0; c < ncol; c++)
            {
                C[r][c] -= min_in_row;
            }

        }
        step = 2;
    }
    // find a zero (Z) in the resulting matrix. If there is no starred
    // zero in its row or column , start Z. Report for each element in the
    // matrix. Go to step 3
    void step_two(int step)
    {
        for (int r =0;r < nrow; r++)
        {
            for (int c = 0; c < ncol; c++)
            {
                if(C[r][c] == 0 && RowCover[r] == 0 && ColCover[c] == 0)
                {
                    M[r][c] = 1;
                    RowCover[r] = 1;
                    ColCover[c] = 1;
                }

            }

        }
        for (int r =0;r < nrow; r++) {
            RowCover[r] = 0;
        }
        for (int c = 0; c < ncol; c++) {
            ColCover[c] = 0;
        }
        step =  3;
    }
    // Cover each column containing a starred zero. If K columns are covered,
    // the starred zeros describe a complete set of unique assignments, In this
    // case, Go to DONE, oterwisem Go to step 4.
    void step_three(int step)
    {
        int colcount = 0;
        for (int r =0;r < nrow; r++) {
            for (int c = 0; c < ncol; c++) {
                if(M[r][c] == 1)
                {
                    ColCover[c] = 1;
                }

            }
            colcount = 0;
        }
        for (int c = 0; c < ncol; c++) {
            if(ColCover[c] == 1)
            {
                colcount += 1;
            }

        }
        if(colcount >= ncol || colcount >= nrow){
            step = 7;
        }
        else {
            step = 4;
        }

    }
    // find a noncovered zero and prime it. If there is no starred zero
    // in the row containing this primed zero, Go to step 5. Otherwise,
    // Continue in this manner until there are no uncovered zeros left.
    // Save the smallest uncovered value and Go to Step 6.
    void step_four(int step)
    {
        int row = -1;
        int col = -1;
        bool done;

        done = false;
        while(!done)
        {
            find_a_zero(row, col);
            if(row == -1)
            {
                done = true;
                step = 6;
            }
            else {
                M[row][col] = 2;
                if( star_in_row(row))
                {
                    find_star_in_row(row, col);
                    RowCover[row] = 1;
                    ColCover[col] = 0;
                }
                else {
                    done = true;
                    step = 5;
                    path_row_0 = row;
                    path_col_0 = col;
                }

            }
        }

    }
    //Construct a series of alternating primed and starred zeros as follows
    //Let Z0 represnt the uncovered primed zero found in step 4. Let Z1 denote
    //the starred zero in the column of Z0 (if any). Let Z2 denote the primed zero
    //terminates at a primed zero that has no stared zero in its column.
    //Unstar each starred zero of the series, star each peimed zero of the series,
    //erase all primes and uncover every line in the matrix. return to step 3
    void step_five(int step)
    {
        bool done;
        int r =-1;
        int c = -1;

        path_count = 1;
        path[path_count - 1][0] = path_row_0;
        path[path_count - 1][1] = path_col_0;
        done = false;
        while (!done)
        {
            find_star_in_col(path[path_count - 1][1], r);
            if(r > -1)
            {
                path_count += 1;
                path[path_count - 1][0] = r;
                path[path_count - 1][1] = path[path_count - 2][1];

            }
            else {
                done = true;
            }
            if(!done)
            {
                find_prime_in_row(path[path_count - 1][0], c);
                path_count += 1;
                path[path_count - 1][0] = path[path_count - 2][0];
                path[path_count - 1][1] = c;

            }
        }
        augment_path();
        clear_covers();
        erase_primes();
        step = 3;
    }
    //Add the value found in step 4 to every element of each coverd row, and subtract
    //it from every elements of each uncovered column. Return to step 4 without
    //altering any stars, primes, or covered lines
    void step_six(int step)
    {
        int minval = MaxValue;
        find_smallest(minval);
        for (int r = 0; r < nrow; r++) {
            for (int c = 0; c < ncol; c++) {
                if(RowCover[r] == 1)
                {
                    C[r][c] += minval;
                }
                if(ColCover[r] == 1)
                {
                    C[r][c] -= minval;
                }
            }
        }
        step = 4;
    }
    void step_seven(int step);
    // methods to support step 4
    void find_a_zero(int row, int col)
    {
        int r = 0;
        int c;
        bool done;
        row = -1;
        col = -1;
        done = false;
        while(!done)
        {
            c = 0;
            while(true)
            {
                if( C[r][c] == 0 && RowCover[r] == 0 && ColCover[c] == 0)
                {
                    row = r;
                    col = c;
                    done = true;
                }
                c += 1;
                if (c >= ncol || done)
                {
                    break;
                }
            }
            r += 1;
            if(r >= nrow)
            {
                done = true;
            }
        }
    }

    bool star_in_row(int row)
    {
        bool tmp = false;
        for (int c =0; c < ncol; c++) {
            if(M[row][c] == 1)
            {
                tmp = true;
            }
        }
        return tmp;
    }

    void find_star_in_row(int row, int col)
    {
        col = -1;
        for (int c = 0; c < ncol; c++) {
            if(M[row][c] == 1)
            {
                col = c;
            }
        }
    }

    // methods to support step 5
    void find_star_in_col(int c, int r)
    {
        r = -1;
        for (int i = 0; i < nrow; i++) {
            if(M[i][c] == 1)
            {
                r = i;
            }
        }
    }

    void find_prime_in_row(int c, int r)
    {
        for (int j = 0; j < nrow; j++) {
            if(M[r][j] == 1)
            {
                c = j;
            }
        }
    }

    void augment_path()
    {
        for (int p = 0; p < path_count; p++) {
            if(M[path[p][0]][path[p][1]] == 1)
            {
                M[path[p][0]][path[p][1]] = 0;
            }else {
                M[path[p][0]][path[p][1]] = 1;
            }
        }
    }

    void clear_covers()
    {
        for (int r = 0; r < nrow; r++) {
            RowCover[r] = 0;
        }
        for (int c = 0; c < ncol; c++) {
            ColCover[c] = 0;
        }
    }

    void erase_primes()
    {
        for (int r = 0; r < nrow; r++) {
            for (int c = 0; c < ncol; c++) {
                if(M[r][c] == 2)
                    M[r][c] = 0;
            }
        }
    }

    // methods to support step 6
    void find_smallest(int minval)
    {
        for(int r = 0; r < nrow; r++)
            for(int c = 0; c < ncol; c++)
                if(RowCover[r] == 0 && ColCover[c] == 0)
                    if(minval < C[r][c])
                        minval = C[r][c];
    }

};

class Hungray:public HungrayAlgorithmBase
{
public:
    Hungray(){}
    ~Hungray(){}

public:
    int n;          //  元素个数
    int *assign;    //  分配结果
    int **mat;      //  代价矩阵
    int **matRcd;   //  代价矩阵
    int totalCost;  //  总成本

public:
    /* 数据读取*/
    bool read(const char* path)
    {
        ifstream extractData;
        extractData.open(path);
        if(!extractData)
            return false;
        extractData>>n;
        assign=new int[n];for(int i=0;i<n;i++)assign[i]=-1;
        mat=new int*[n];for(int i=0;i<n;i++)mat[i]=new int[n];
        matRcd=new int*[n];for(int i=0;i<n;i++)matRcd[i]=new int[n];
        for(int i=0;i<n;i++){
            for(int j=0;j<n;j++){
                extractData>>mat[i][j];
                matRcd[i][j]=mat[i][j];
            }
        }
        extractData.close();
        totalCost=0;
        return true;
    }

    /* 行归约 */
    void rowSub()
    {
        int *minEmt=new int[n];for(int i=0;i<n;i++)minEmt[i]=(int)1e8;

        for(int i=0;i<n;i++)for(int j=0;j<n;j++)if(mat[i][j]<minEmt[i])minEmt[i]=mat[i][j];

        for(int i=0;i<n;i++)for(int j=0;j<n;j++)mat[i][j]-=minEmt[i];

        delete []minEmt;
    }

    /* 列归约 */
    void columnSub()
    {
        int *minEmt=new int[n];for(int j=0;j<n;j++)minEmt[j]=(int)1e8;

        for(int j=0;j<n;j++)for(int i=0;i<n;i++)if(mat[i][j]<minEmt[j])minEmt[j]=mat[i][j];

        for(int j=0;j<n;j++)for(int i=0;i<n;i++)mat[i][j]-=minEmt[j];

        delete []minEmt;
    }

    /* 检验最优 */
    bool isOptimal(int *assign)
    {
        int *tAssign=new int[n];
        for(int i=0;i<n;i++)
            tAssign[i]=-1;
        int *nZero=new int[n];
        bool *rowIsUsed=new bool[n];
        bool *columnIsUsed=new bool[n];

        for(int i=0;i<n;i++)
            rowIsUsed[i]=columnIsUsed[i]=0;
        int nLine=0;
        while(nLine<n){
            for(int i=0;i<n;i++)
                nZero[i]=0;
            for(int i=0;i<n;i++){
                if(rowIsUsed[i]==1)
                    continue;
                for(int j=0;j<n;j++){
                    if(columnIsUsed[j]!=1&&mat[i][j]==0)
                        nZero[i]++;
                }
            }

            int minZeros=n;
            int rowId=-1;
            for(int i=0;i<n;i++){
                if(rowIsUsed[i]==0&&nZero[i]<minZeros&&nZero[i]>0){
                    minZeros=nZero[i];
                    rowId=i;
                }
            }
            if(rowId==-1)
                break;
            for(int j=0;j<n;j++){
                if(mat[rowId][j]==0&&columnIsUsed[j]==0){
                    rowIsUsed[rowId]=1;
                    columnIsUsed[j]=1;
                    tAssign[rowId]=j;
                    break;
                }
            }
            nLine++;
        }
        for(int i=0;i<n;i++)
            assign[i]=tAssign[i];


        delete []tAssign;
        delete []nZero;
        delete []rowIsUsed;
        delete []columnIsUsed;

        for(int i=0;i<n;i++)if(assign[i]==-1)
            return false;

        return true;
    }

    /* 矩阵变换 */
    void matTrans()
    {
        bool *rowTip=new bool[n];
        bool *columnTip=new bool[n];
        bool *rowLine=new bool[n];
        bool *columnLine=new bool[n];

        for(int i=0;i<n;i++)
            rowTip[i]=columnTip[i]=rowLine[i]=columnLine[i]=0;
        //打勾
        for(int i=0;i<n;i++)
            if(assign[i]==-1)
                rowTip[i]=1;

        while(1){
            int preTip=0;
            for(int i=0;i<n;i++)
                preTip=preTip+rowTip[i]+columnTip[i];
            for(int i=0;i<n;i++){
                if(rowTip[i]==1){
                    for(int j=0;j<n;j++){
                        if(mat[i][j]==0)columnTip[j]=1;
                    }
                }
            }
            for(int j=0;j<n;j++){
                if(columnTip[j]==1){
                    for(int i=0;i<n;i++){
                        if(assign[i]==j)rowTip[i]=1;
                    }
                }
            }
            int curTip=0;
            for(int i=0;i<n;i++)
                curTip=curTip+rowTip[i]+columnTip[i];
            if(preTip==curTip)
                break;
        }

        /* 画线 */
        for(int i=0;i<n;i++){
            if(rowTip[i]==0)
                rowLine[i]=1;
            if(columnTip[i]==1)
                columnLine[i]=1;
        }

        /* 找最小元素 */
        int minElmt=(int)1e8;
        for(int i=0;i<n;i++)
            for(int j=0;j<n;j++)
                if(rowLine[i]==0&&columnLine[j]==0&&mat[i][j]<minElmt)
                    minElmt=mat[i][j];
        /* 变换 */
        for(int i=0;i<n;i++)
            if(rowTip[i]==1)
                for(int j=0;j<n;j++)
                    mat[i][j]-=minElmt;
        for(int j=0;j<n;j++)
            if(columnTip[j]==1)
                for(int i=0;i<n;i++)
                    mat[i][j]+=minElmt;

        delete []rowTip;
        delete []columnTip;
        delete []rowLine;
        delete []columnLine;
    }

    /* 匈牙利算法 */
    void runMunkres()
    {
        read("/home/freja/lxy/GITEE/data-structures-and-algorithms/DataStructuresAndAlgorithms/Hungary/cost.dat");//读取数据
        rowSub();//行归约
        columnSub();//列归约
        //如果不能找到n个独立的0元素，则对矩阵进行变换
        while(!isOptimal(assign)){
            matTrans();
        }

        for(int i=0;i<n;i++)
            totalCost+=matRcd[i][assign[i]];
        for(int i=0;i<n;i++)
            delete []mat[i];
        delete []mat;
        for(int i=0;i<n;i++)
            delete []matRcd[i];
        delete []matRcd;
    }
};


int testMunkre()
{
    Hungray *p;
    p = new Hungray;
    // 调用匈牙利算法
    p->runMunkres();
    // 输出结果
    cout<<"总成本为"<<p->totalCost<<endl;
    for(int i=0;i<p->n;i++)
        cout<<"员工"<<i+1<<"-->"<<"任务"<<p->assign[i]+1<<endl;

    cin.get();
}


#endif // HUNGRAY_H
